Boundary conditions for augmented plane wave methods 
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The augmented plane wave method uses the Rayleigh-Ritz principle for basis functions that 
are continuous but with discontinuous derivatives and the kinetic energy is written as a pair of 
gradients rather than as a Laplacian. It is shown here that this procedure is fully justified from the 
mathematical point of view. The domain of the self-adjoint Hamiltonian, which does not contain 
functions with discontinuous derivatives, is extended to its form domain, which contains them, and 
this modifies the form of the kinetic energy. Moreover, it is argued that discontinuous basis functions 
should be avoided. 
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I. INTRODUCTION 

The augmented plane wave (APW) method and its by- 
products are among the most powerful electronic struc- 
ture calculation schemes 1 -. The augmented plane wave 
method was invented by Slater in 1937—. Its basic idea is 
to partition the crystal into two types of regions: (i) non- 
overlapping spheres S a centered around all constituent 
atomic sites and (ii) the remaining interstitial region. 
The solutions ip( r ) of the Schrodinger equation are ex- 
panded over basis wavefunctions fa(r) which are plane 
waves in the interstitial region and spherical wavea^ in- 
side each sphere. The coefficients of this expansion are 
obtained through the Rayleigh-Ritz principle 4 ^ by find- 
ing the extrema of the Raylcigh quotient— 
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(1) 



A subtle point of the APW method is the boundary 
conditions for the basis functions fa at the surface of the 
spheres S a . Various choices can be found in the litera- 
ture. In the linearized APW method^, the basis func- 
tions and their derivatives are continuous at the spheres. 
In the original APW method, as well as in the recent 
linearized APW method with local orbitals^, the basis 
functions are only continuous at the spheres. Moreover, 
many authors proposed to use basis functions that are 
discontinuous at the spheresiSiii. 

The boundary conditions are important to define the 
basis functions but also to define the Rayleigh quotient 
itself. If the functions and their derivatives are contin- 
uous, an integration by parts (or Green's first identity) 
shows that the kinetic energy term of the Hamiltonian 
can be written in two equivalent forms: 



i*(r)A^-(r)dr 



(V&(r))*.V^-(r)dr. (2) 



If the derivatives of the basis functions are not contin- 
uous, then the left-hand side (i.e. the Laplacian form) 



and the right-hand side (the gradient form) of equation 
P)l are generally different and the question is: which 
one (if any) is correct? Slater's answer is "the second, 
more fundamental form is the correct one, as it is the one 
which directly enters the variation principle from which 
Schrodinger's equation is derived"—. Obviously, this is 
not a very convincing argument because the variational 
principle can be written with the kinetic energy in the 
Laplacian form as well as in the gradient form. Still, this 
argument (and practical success) has been up to now the 
only one to justify the use of the gradient form. If the 
basis functions are not continuous, additional terms are 
introduced 1011 - but the method looses many of the at- 
tractive features of the standard APW approach. 

The purpose of this paper is to show that the question 
of the boundary conditions in the APW method is not 
a matter of physical principles, but a well defined math- 
ematical question with a precise answer which is: basis 
functions with discontinuous derivatives are valid if the 
gradient form is used and discontinuous basis functions 
should be employed with much caution because they do 
not allow for a correct variational principle. 

To obtain this answer, we make a short trip in the 
mathematical theory of Hamiltonian operators. We first 
consider the fact that realistic Hamiltonians H are de- 
fined on a subspace D(H) of the Hilbert space called the 
domain of H. The Rayleigh-Ritz principle holds for basis 
functions in the domain of H. However, the domain of 
H does not contain functions with discontinuous deriva- 
tives. To extend the possible basis functions, we use the 
form domain Q(H) of H. We show that the Rayleigh- 
Ritz principle is still valid on Q(H) if the gradient form 
of the kinetic energy is used. These results fully confirm 
Slater's intuitions. Finally, we show that the Rayleigh- 
Ritz principle is not valid if the basis functions are not 
continuous. This implies a number of unpleasant con- 
sequences which show that discontinuous basis functions 
should be used with great care. 
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II. THE HAMILTONIAN AND ITS DOMAIN 

It is sometimes thought that a Hamiltonian is a Her- 
mitian operator H defined on a Hilbert space H in the 
sense that for any (j>, -0 in H we have (H(f>\ijj) = (<f)\H%li). 
This not true for realistic Hamiltonians. The Hellinger- 
Toeplitz theorem (ref. P-84) tells us that an operator 
H can satisfy (H<fi\ip) = (cf)\Hip) for all 0, ip in H if 
and only if the eigenenergies of H are between a mini- 
mum and a maximum finite energy. Realistic Hamitoni- 
ans have a minimum energy (that of the ground state) 
but they have no maximum energy. From the physical 
point of view, a system represented by a Hamiltonian 
with a maximum energy would completely reflect pho- 
tons or electrons above that energy. Experiment tells us 
that very high energy photons and electrons go through 
matter so such Hamiltonian is not realistic. From the 
mathematical point of view, the Hamiltonians made of 
a kinetic energy term plus a potential energy term that 
is periodic or go to zero at infinity (plus some regularity 
assumptions) has no maximum energy. 

As a consequence, a realistic Hamiltonian H is not de- 
fined on a Hilbert space Ti but on a subspace D(H) of 
Ti., called the domain of H, such that H \(j>) € Tt for any 
4> 6 D(H). We shall see that the question of the bound- 
ary conditions turns out to be a question concerning the 
domain of H . 



A. Self-adjoint Hamiltonian 

Since the difference between the Laplacian and the gra- 
dient forms comes from the kinetic energy, we shall take 
the example of H = — A in the Hilbert space of square in- 
tegrable functions L 2 (R 3 ). What is the domain of HI In 
fact, many domains are possible. Let us be very conserva- 
tive and take the domain Z3 m i n to be the set of infinitely 
differentiable functions with compact support (i.e. which 
are zero outside a bounded domain of R 3 ). It is clear that 
£>min in included in L 2 (R 3 ) and that H is well defined 
on -D m i n . However this is not the domain that we want 
because the domain of the adjoint of H is not -D m in- 

The Hamiltonian H with domain D min has the ad- 
joint -fft with domain D(H^) defined as follows (ref. 
El p.252). For any V € D(W) and <f) G Anin, 
/ f(r)(5V(r))dr = - / (A0(r))V(r)dr. The func- 
tion (j) is smooth and with compact support, so the pre- 
vious relation shows that H^ip(r) = — A-0(r) in the sense 
of distributions. As such, ip does not need to be contin- 
uous for Atp(r) to have a meaning in the sense of distri- 
butions. But W is an operator from H to H, therefore 
Aip(r) must belong to Ti. To summarize, the adjoint of 
H with domain D m i n is the operator —A acting in the 
sense of distributions, and its domain is the set of square 
integrable functions ip such that Aip is square integrable. 
It can be shown that the domain of is the second 
Sobolev space! 3 . W 2 (ref. Q p. 55). 

From this example, it is clear that the domain of H and 



of can be different. An operator is called self-adjoint 
when H = H* and the domain of H is equal to the 
domain of H* . We want self-adjoint operators because 
many properties that we take for granted in quantum 
physics (Rayleigh-Ritz principle, expansion over eigen- 
states, unitarity of exp(itH), etc.) are valid only for self- 
adjoint operators (ref. [H p. 256). 

Is there a domain over which H — —A is self-adjoint? 
Yes, and this domain is unique, it is W2 (ref. [lj, p. 54). 
This unicity is not general. For example, the radial 
Schrodinger operator 

h = d " - - I (3) 
dr 2 r dr r 2 ' 

has only one self-adjoint domain when t ^ 0, but has a 
family of self-adjoint domains when I = 0i£ (technically, 
this is because the Laplacian is not essentially self-adjoint 
in a space of dimension n < 4, see~ p. 129). 



B. Rayleigh-Ritz principle 

The importance of the self-adjoint domain D(H) of 
H for the APW method comes from the fact that, if 
H is self-adjoint and bounded from below (i.e. has a 
ground state) , then for any n-dimensional subspace V of 
D(H) with orthonormal basis 1 - <f>i, the fc-th eigenvalue 
Afe of the matrix (4>i\H\<fij) is larger than the fc-th eigen- 
value Ek of H. Moreover, for realistic Hamiltonians, it 
can be shown that Ao converges to the ground state en- 
ergy Eq of H as the dimension n of V increases (ref. ^3 
p. 82). In favorable cases all the eigenvalues Afc converge 
to Ek- This happens in particular when the eigenvalues 
are discret o 18 ! 19 , which is the case of the Bloch energies 
E n (\t) for a given Brillouin zone vector k. 

In other words, if the basis functions <fii belong to 
D(H), the Rayleigh-Ritz principle gives an approximate 
ground state with an energy greater than the true ground 
state, and the approximate ground state converges to the 
true ground state by increasing the number of basis func- 
tions. This is exactly the properties we assume when us- 
ing the Rayleigh-Ritz principle. It seems that the ques- 
tion of the boundary conditions is solved by determining 
which basis functions belong to D(H). However, if the 
derivative of a basis function (f> is not continuous at the 
sphere boundary, then the second derivative involved in 
the condition A(j> £ L 2 (R 3 ) brings a Dirac S function at 
the sphere radius, and the square of the S function is not 
defined, so A^ does not belong to L 2 (M 3 ). Therefore, for 
a basis function to be in D(H), it must be continuous 
with continuous derivatives at the sphere boundary. 



III. RELAXING THE BOUNDARY CONDITION 

The previous result is disappointing because the con- 
dition of a continuous derivative at the boundary seems 
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much too strong. Many successful calculations of solid- 
state properties were carried out with basis functions 
that are only continuous at the sphere boundary. There 
should be a way to enlarge the domain of the basis func- 
tion. This w_ay exists, it is called the theory of quadratic 
forms (ref. Ojp.276). 

A. Quadratic forms 

We consider only the quadratic form associated with 
a self-adjoint operator H, and more precisely with H = 
—A. This quadratic form is defined by 

Q(6l/0 = (tt>\H^) = - /V(r)AV>(r), (4) 

for any (f> and tp in D(H). In general, Q can be extended 
to a domain Q(H) larger than D(H). To increase the 
domain of Q we rewrite the quadratic form as 

W» = / V<f (r)-W(r). (5) 

The transformation from Q to J5J is valid if <fi and ip 
are in D(H). Now the form j£| has a larger domain 
Q(H) defined as the set of 0(r) e L 2 (R 3 ) such that 
V</>(r) 6 L 2 (R 3 ). This is a Hilbert space called the first 
Sobolev space W\ (ref. LLJ p. 50). This way of enlarging 
the domain of Q is mathematically unambiguous: equa- 
tion (0J defines a closable quadratic form on D(H) and 
its (unique) closure is the quadratic form (3 given by 
equation © with domain Q{H) = Wi (ref. HI d. 177). 

B. Rayleigh-Ritz principle for quadratic forms 

This extension of the domain of H is interesting only if 
we can use the Rayleigh-Ritz principle. It is not difficult 
to show that this is the case 2 ^: For any n-dimensional 
subspace V of Q(H) with orthonormal basis <fii, the k-th 
eigenvalue of the matrix (<j)i\H\<pj) is larger than the 
A;-th eigenvalue Ek of H. In favorable cases, Afc converge 
to Ek by increasing the dimension of V—. 

We can now determine the regularity of the basis func- 
tions. If they are continuous at the sphere boundary, the 
derivative along the sphere radius has a step, which is 
a square integrable function. Therefore, the APW ba- 
sis functions are in Q(H) if they are continuous. We 
can even prove that they have to be continuous. This 
is done with Sobolev's lemma (ref. ^] p. 52). How- 
ever, we must not consider the Sobolev space VFi(R 3 ) 
but the Sobolev space Wi(R) because the use of radial 
wave functions transforms the three-dimensional gradi- 
ent into a one-dimensional gradient. Inside sphere S a , 
the radial function is /if m (r Q ) (see note0), outside this 
sphere an expansion of the planewave over spherical har- 
monics shows that the radial function is proportional to 
the Bessel function ji(kr a ). We want the gradient of 



this radial function to be square integrable, so the ra- 
dial function is in W\ (R). Sobolev's lemma tells us now 
that all functions of Wi(R) are continuous, so the radial 
functions are continuous in Q(H). 



C. Realistic Hamiltonians 

In the foregoing sections, we dealt only with the op- 
erator —A, whereas interesting Hamiltonians have the 
form H = —A + V, where V is a real potential, which 
is smooth except for Coulomb singularities at the nuclei. 
In an atom or a molecule, V goes to zero at infinity, in a 
crystal V is periodic. It turns out that the addition of a 
potential does not modify our conclusions. 

For an atom or a molecule, the potential V is the sum 
of a finite number of Coulomb potentials plus a bounded 
and square integrable potential (e.g. the exchange and 
correlation potentials in the local density approximation 
to the density functional theory). For these potentials, 
H is self-adjoint on the same domain as —A (ref. [lj 
p. 165). Thus, the Coulomb singularities do not change 
the domainSi. Moreover, the form domain is the same 
as that of —A (this is essentially because the difference 
between D(H) and Q(H) comes from the fact that the 
first one uses A<f> whereas the second one uses V</>, and 
this does not depend on V if the singularities of V are not 
stronger than the Coulomb potential) . Therefore, the 
whole discussion is still valid for realistic Hamiltonians 
representing an atom or a molecule. 

For a solid, similar results are available. If £1 is a 
primitive cell of the crystal, the Hilbert space 7ik for 
the Brillouin zone vector k is the set of functions /(r) 
that are square integrable on and satisfy the Bloch 
conditions /(r + R) = exp(ik • R)/(r) and V/(r + R) = 
exp(zk • R) V/(r), for any point r on the surface of fl and 
any lattice vector R such that r + R is also on the surface 
of fi (ref. [l]Jp. 303). If HQs) is the Hamiltonian ~A + V 
acting on and if the periodic potential V is smooth 
except for a finite number of Coulombic singularities, the 
domain of H(k) is the Sobolev space W2(0), i.e. the set 
of functions <j> of Hk such that Atfi is square integrable on 
O. The Bloch conditions ensure that the quadratic form 
associated with H (k) is 

W» = - / r(r)AV(r)+ f 0*(r)l/(r)V(r) 

= { V0*(r) ■ W(r) + / 4>*(r)V(r)iP{r). 
Jn Jn 

for any <f> and ip in D(H(k)). Again the second ex- 
pression for Q can be extended to the form domain 
Q(H(k)) = Wi(Q). Therefore, the whole discussion re- 
mains valid in the presence of a periodic potential and the 
basis functions must be continuous at the sphere bound- 
ary. 
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D. Discontinuous functions 

We can try to further relax the boundary conditions 
and to take basis functions which are not continuous at 
the sphere boundary. This was proposed by a number of 

]u ^ : n 10.11.23.24.25.26.27.28.29.30.31.32.33.34.35.36.37.38.39 

Among practitioners of APW methods, the prevailing 
opinion is that discontinuous boundary conditions are le- 
gitimate. For example Shaughnessy and coll— state that 
"There is no a priori reason why the different forms of 
the basis functions, inside and outside the spheres, should 
match in either value or slope at the sphere surface" . 

We saw that there is a good reason for using continuous 
basis functions. This ensures the validity of the Rayleigh- 
Ritz principle, from which we can deduce that to each 
state of the eigenvalue problem on the reduced vector 
space V corresponds a state of H, that the eigenvalues of 
the reduced problem are larger than the corresponding 
eigenvalue of H and that convergence to the eigenvalues 
of H can be achieved by increasing the dimension of V. 

One could be optimistic and hope that the Rayleigh- 
Ritz principle is still valid for discontinuous functions. 
However, such a hope is generally not fulfilled. This was 
shown by Brownsteinifi who investigated simple Hamilto- 
nians H with discontinuous basis functions. He showed 
that discontinuity can lead to rather unpleasant conse- 
quences: (i) spurious states appear, for example states 
with energies smaller than the ground state of H, (ii) 
increasing the dimension of the reduced vector space V 
does not always improve the agreement with the true 
eigenvalues, (iii) the degeneracy of the states can be 
wrong. Moreover, when complex wavefunctions are used 
(this is necessary in band-structure calculations), some 
formulations give an energy which is a complex number—. 
These phenomena are obviously quite dangerous when 
one tries to calculate (or predict) properties of matter 
and it seems reasonable that discontinuous basis func- 
tions are not used in practical implementations of the 
APW method. Of course, as shown by Brownstein, cal- 
culations using discontinuous basis functions can yield 
correct results in favorable cases. But no general princi- 
ple ensures that the results of the computation are the 
correct energies and eigenstates of H. 

IV. CONCLUSION 

The results of this paper are clear-cut: for the 
Rayleigh-Ritz principle to be valid the basis functions 
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